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This paper presents selected results from an ongoing effort to develop an aerodynamic database from 
Reynolds-Averaged Navier-Stokes (RANS) computational analysis of airfoils and wings at stall and 
post-stall angles of attack. The data obtained from this effort will be used for validation and re- 
finement of a low-order post-stall prediction method developed at NCSU, and to fill existing gaps in 
high angle of attack data in the literature. Such data could have potential applications in post-stall 
flight dynamics, helicopter aerodynamics and wind turbine aerodynamics. An overview of the NASA 
Te trUSS CFD package used for the RANS computational approach is presented. Detailed results for 
three airfoils are presented to compare their stall and post-stall behavior. The results for finite wings 
at stall and post-stall conditions focus on the effects of taper-ratio and sweep angle, with particular at- 
tention to whether the sectional flows can be approximated using two-dimensional flow over a stalled 
airfoil. While this approximation seems reasonable for unswept wings even at post-stall conditions, 
significant spanwise flow on stalled swept wings preclude the use of two-dimensional data to model 
sectional flows on swept wings. Thus, further effort is needed in low-order aerodynamic modeling of 
swept wings at stalled conditions. 


Nomenclature 

AR 

= 

aspect ratio 

a 

= 

VGRID stretching/growth factor 

b 

= 

VGRID stretching/growth factor 

CFD 

= 

computational fluid dynamics 

c d 

= 

airfoil or section drag coefficient 

Cf 

= 

skin friction coefficient 

C L 

= 

wing lift coefficient 

Ci 

= 

airfoil or section lift coefficient 

r 

K -'m 

= 

airfoil or section pitching moment coefficient about the quarter-chord 

c P 

= 

pressure coefficient 

C ref 

= 

reference chord length 

Re 

= 

Reynolds number based on airfoil chord 

SA 

= 

Spalart-Allmaras turbulence model 

TE 

= 

trailing edge 

At 

= 

physical time step 

Uoo 

= 

Freestream reference velocity 

y+ 

= 

non-dimensional wall distance 

A zi 

= 

height of the I th viscous layer 

a 

= 

angle of attack, deg 

A 

= 

angle of sweepback, deg 

A 

= 

taper ratio 
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I. Introduction 


Airfoil lift and moment data in the low angle of attack regime is readily available from a multitude of sources- from 
Abbott and von Doenhoff [1] and those from the University of Stuttgart [2] and University of Illinois at Urbana- 
Champaign [3,4,5], to modern computational approaches designed to predict sectional aerodynamic characteristics 
based on arbitrary input geometry, such as XFOIL [6]. For many applications, data in this linear regime is sufficient. 
However, fields such as wind turbine aerodynamics, helicopter aerodynamics and post-stall flight dynamics of fixed- 
wing aircraft require data to extend beyond aerodynamic stall. Efforts have been made in the wind turbine community 
and the helicopter aerodynamics community to extend airfoil data into the post-stall regime. In both wind turbine and 
helicopter aerodynamics, local blade sections close to the root of the rotor may experience very high angles of attack. 
Models for airfoil force and moment coefficients at high angle of attack conditions have been developed experimentally 
[7,8], which has led to researchers proposing empirical models based on flat plate theory [9]. The empirical models 
developed from experiment require that both the maximum C/ and the corresponding a sta u at which this lift coefficient 
occurs be known reliably before theoretical flat plate data may be fitted to extend the data well into post-stall. 

Very little data exists in the literature covering the stall behavior of finite wings, especially that which extends 
deep into post-stall. Published work in post-stall wing aerodynamics often covers general stall behavior, dependent 
on factors such as planform, without providing detailed force and moment data beyond initial stall [10]. Some studies 
propose empirical methods based on theory to extend existing force and moment coefficient data-sets deep into post- 
stall taking into account 3D effects with some correction for aspect ratio [8]. An interesting, purely experimental study, 
that covers the stall of both a 2D airfoil section and 3D wings of various aspect ratios was performed by Ostowari and 
Naik in 1985 [1 1]. The study presented consistent lift coefficient versus angle of attack data for various NACA 44XX 
series airfoils and 3D rectangular wings with a range of aspect ratios having the same airfoils as cross sections. 

A database of post-stall airfoil and wing data, with a similar scope to the Ostowari and Naik study, was desired 
partly to address this dearth of data in the literature, but primarily to support a local effort within the Applied Aero- 
dynamics Group at North Carolina State University. The ongoing effort involves developing a low-order model of 
post-stall aerodynamics for finite wings via use of existing linear low-order methods (VLM, Weissinger or LLT) cor- 
rected for nonlinear sectional airfoil behavior. Corrections are accomplished via a decambering approach used to 
mimic the nonlinear aerodynamics caused by flow separation. Details of the development of the low-order method 
may be found in Ref. [12], with the current status pertaining to its use in real-time simulation of aircraft flight dynam- 
ics described in Ref. [13]. In addition to geometry information, the low-order post-stall model requires sectional 2D 
airfoil data as input for each of the control points. This sectional input data (C/ - a, Cd~ a, and C m - a curves) defines 
the convergence criteria for the low-order post-stall calculation for finite wings. Outputs include the total aerodynamic 
force and moment coefficients and span wise distributions of these coefficients. 

Figure 1 describes the relationship between the low-order method and the higher order Computational Fluid Dy- 
namics (CFD) effort, and how the results from the latter are to be used for initial validation and further refinement of 
the former. Two dimensional CFD is performed on airfoil sections, and the outputs are used as inputs to the low-order 
method. The 3D geometries utilizing the same 2D cross section are run in both CFD and the low-order method.The 
results are then compared on the basis of total force/moments and span wise force/moment distributions. CFD solu- 
tions and accompanying flow visualizations allow for further development and refinement of the low-order method. 
Comparisons between the low-order method and CFD are beyond the scope of this paper. This paper aims to discuss 
the methodology used in generating CFD solutions of airfoils and wings in post-stall angles of attack and to present 
the results in this regime as dependent on geometry and flow conditions. 

This paper begins by describing the CFD software package that was chosen to generate consistent flow solutions 
for airfoils and wings. The methodology developed to effectively use the CFD package is presented next - from 
generating usable geometry, creating a suitable grid over the computational domain, running the flow-solver, checking 
solution convergence, and post-processing to extract desired quantities from the outputs. The 2D and 3D results are 
presented which show effects from geometry, including wing sweep, and flow parameters. Selected results are shown 
in more detail using flow visualization and/or span wise local lift coefficient distributions. 
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Figure 1: Approach for using CFD for initial validation and further refinement of the low-order method. Blue empha- 
sizes that a consistent grid spacing is used with the 2D and 3D CFD discretizations. 


II. Background: NASA Tetrahedral Unstructured Software System 

The NASA Tetrahedral Unstructured Software System (TetrUSS) [14] is a package of loosely integrated software, 
developed to allow rapid aerodynamic analysis from simple to complex problems. The system has its origins in 1990 at 
the NASA Langley Research Center and has won the NASA Software of the Year Award twice. TetrUSS has been used 
on high priority NASA programs such as Constellation and the new Space Launch System for aerodynamic database 
generation along with work in the Aviation Safety Program. The component software packages are assembled such 
that a user follows a systematic process during each application of TetrUSS. There are software packages for preparing 
geometries for grid generation (GridTool), generating unstructured surface and volume grids (VGRID) and calculating 
flow solutions (USM3D). Post-processing the solutions with TetrUSS can be done using the included Simple View 
software or by easily converting for use with other commercial packages (eg Tecplot, EnSight etc.). 

For preparing geometries for grid generation, GridTool is used to generate the necessary VGRID [15] input files. 
GridTool can read Non-uniform Rational B-Spline (NURBS) curves and surfaces through an Initial Graphics Exchange 
Specification (IGES) file, as well as PLOT3D point cloud surface definitions. The geometric surfaces are then defined 
by way of surface patch construction. Each patch has a specified boundary condition such as a viscous or inviscid 
surface and a family definition for users to group related patches together. Grid spacing parameters are also defined 
and controlled within GridTool. Sources are placed in three dimensional space by a user in order to control the size 
and growth rates of the tetrahedral cells. Numerous classes of sources are available to control the grid topography. 
Nodal sources and line sources are typically used in most cases, while volume sources are available for use in specific 
cases requiring control over a large volume of the domain. Other parameters defined in GridTool are the viscous layer 
spacing parameters and the maximum and minimum tetrahedral sizes. 

VGRID is the unstructured grid generation tool used in the TetrUSS Package. Viscous layer generation is ac- 
complished via the Advancing Layers Method (ALM) [16]. Tetrahedral cells are generated in an orderly manner, 
“marching” nodes away from the surface. The size and growth of these cells is controlled by Equation 1. 
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( 1 ) 


Az i+ i = Azi[l + a(l + b)‘] 

In this equation, the height of the I th layer is determined by an initial spacing parameter, Azu and two stretching/growth 
factors a and b. Once the height of the i th layer reaches the size of the background sources specified by the user in 
GridTool, no more cells are formed and viscous layer generation is complete. After the viscous layers are generated, 
VGRID then utilizes the Advancing Front Method (AFM) [17] for the generation of the inviscid portion of the volume 
grid. VGRID can not always close the grid completely. When this occurs, a slower but more robust auxiliary code 
called POSTGRID is used to complete the formation of the remaining tetrahedral cells. 

The flow solver at the core of the TetrUSS package is USM3D [18]. USM3D is a parallelized, tetrahedral cell- 
centered, finite volume Reynolds Averaged Navier- Stokes (RANS) flow solver. It computes the finite volume solution 
at the centroid of each tetrahedral cell and utilizes several upwind schemes to compute inviscid flux quantities across 
tetrahedral faces. USM3D has numerous turbulence models implemented for use; the Spalart-Allmaras (SA) one- 
equation model and Menter Shear Stress Transport (SST) two equation model were used in this study. Some additional 
capabilities that USM3D has implemented are dynamic grid motion and overset grids. 

III. Methodology: Developing the Aerodynamic Database 

The aerodynamic database is desired to have high fidelity flow solutions for a wide variety of 2D airfoils and 3D 
geometries. Flow solutions would include a large range of angles of attack to encompass pre-stall, stall and post- stall 
flow regimes. The data gained from these simulations will be vital in assisting the further development and validation 
of the low-order method mentioned in the introduction as well as to fill gaps in the currently available high angle of 
attack aerodynamic data for arbitrary geometries. An efficient process to go from a geometry to a converged flow 
solution was developed and is discussed in the following sub-sections. 

A. Geometry Generation 

Traditional Computer Aided Design (CAD) software would be more than adequate for the creation of the desired 
geometries, however these tools are not geared specifically towards the modeling of wings and airfoils. Understanding 
this, the recently released parametric modeling tool, Open Vehicle Sketch Pad (OpenVSP) [19], was chosen as the 
geometry generation tool. A flow chart showing the process of geometry generation can be seen in Figure 2. OpenVSP 
is a modeling package developed and released by NASA Langley Research Center in Hampton, Virginia. The unique 
concept that OpenVSP provides is that it allows a user to drag and drop generic aircraft components (such as a wing) 
into the modeling area, and directly manipulate familiar geometric parameters. Consequently, it is simple to insert a 
wing, change its root chord, tip chord, span, etc. and view the resulting geometry in real time. Aerodynamic reference 
quantities can also be automatically calculated for the user. Airfoils cross section generation is also simplified. A 
user can select any 4 or 5 digit series NACA airfoil or load in a formatted airfoil coordinate file for use on any lifting 
surface. 

The 3D wing and corresponding airfoil for each case to be analyzed were generated in OpenVSP. In order to read 
the geometry into GridTool, the file must be in the IGES format. Vehicle Sketch Pad does not output IGES files, 
thus each geometry must be exported as a Rhino3D formatted file. The Rhinoceros NURBS modeling package was 
used to convert the geometry into the necessary IGES file as well as make small modifications to the geometry. Some 
grid generation failures were encountered due to how OpenVSP closes the trailing edge of the wing/airfoil geometries 
(it always forces a sharp trailing edge). In some cases the sharpness had to be removed to ensure successful grid 
generation. 
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Vehicle Sketch Pad to 
generate CAD surface 


OpenWE 



Rhino 3D for conversion 
to IGES format 




Generate grid with 
GridTool/VGRID 


Rhinoceros® 



NURBS modeling for Windows 



Figure 2: Process of geometry generation to grid generation 


B. Grid Generation 

Grid generation parameters were generalized such that, between different geometries, parameters such as source place- 
ment and viscous spacing had minimal required changes. This meant that from initial geometry generation to a com- 
pleted grid would require only a matter of hours. Establishing this commonality and routine for grid generation enabled 
the generation of adequate grids for many configurations in a short time span. An example placement of sources for 
a simple tapered wing is shown in Figure 3. A series of line sources are utilized, in the span wise direction of each 
wing at differing chord- wise locations. Anisotropic stretching [20] as high as 10:1 was used near the root of the wing, 
transitioning to isotropic cells near the wing tips. For viscous tetrahedral layer generation, the height of the first layer 
(Azi) is Reynolds number dependent. A viscous layers spacing tool called USGUTIL, was used to determine the 
height of the initial viscous layer. In order to have adequate number of cells in the viscous layers, the y + of the first 
node was set to be 3, this would ensure that the y + of the first cell center would be less than 1 (approx. 0.75) as is 
required for a fully viscous Navier-Stokes solution. The values used for the grid growth parameters (« a and b) in Eq. 
1 were 0.15 and 0.02 respectively [20]. A grid sensitivity study was performed on a rectangular wing with aspect 
ratio of 12 to determine adequate grid sizing. It was found that a grid sizing of 5-9 million tetrahedral cells showed 
changes in Cl, max of approximately 0.02 between the grids. This method of grid generation was applied to all 3D wing 
geometries with the typical grids averaging between 9-12 million tetrahedral cells. 

For airfoil calculations, a quasi 2D grid was generated on a constant-chord, short-span wing, between two reflection 
plane boundary condition patches. Figure 4 shows a completed grid for an NACA 0012 airfoil. A general goal was 
set to maintain very similar grid density between the 3D wing grids and the airfoil grids. This is necessary because 
the airfoil results were being used as input data into the low-order method discussed in the introduction while the 3D 
wing results from USM3D were being used to assess the accuracy of the low-order method (Figure 1). Therefore a 
separate grid sensitivity study was not performed specifically for airfoils. Typical airfoil grids were on the order of 
300,000 tetrahedral cells and were generated using a nearly identical source placement as the 3D wing. 

C. Flow Solution Generation 

All solutions with the USM3D solver were computed with time-accurate Reynolds-Averaged Navier-Stokes (RANS). 
The limitations of RANS for modeling massively separated flows are well known. The more preferred Detached Eddy 
Simulation (DES) modeling will provide better physical representation of 3D separated flow, but with an order-of- 
magnitude more expense. Since this investigation requires generation of many flow solutions, the initial focus is to 
determine if time-accurate RANS can provide sufficient engineering accuracy for capturing the salient aerodynamic 
characteristics of wings at stall and post-stall conditions. Furthermore, a consistent modeling is desired between the 
2D airfoils and 3D wings. 

All computations were advanced at a characteristic time step of AC = At • Uoolc re f = 0.02 using a second or- 
der time-accurate scheme with three-point backward differencing and physical time stepping. The number of sub- 
iterations for each time step was set to between 10 and 15 to ensure adequate sub-iteration convergence. The Spalart- 
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Figure 3: Screen capture of GridTool source placement 
on a tapered wing. 



Figure 4: Screen capture of airfoil grid. 


Allmaras (SA) [21] one equation turbulence model was used almost exclusively, however some simulations were 
performed with the two equation SST turbulence model to understand the difference in final solution quantities. The 
solver was run on both a NASA Langley computer cluster and the North Carolina State University High Performance 
Computing (NCSU HPC) cluster. Making use of USM3D’s parallel computation capabilities, each grid was parti- 
tioned into 28-64 equal zones which could be loaded onto 28-64 individual processors, reducing calculation times 
significantly. To further increase productivity, a series of Unix scripts were developed to generate the required input 
files for job submission and minor post-processing of completed jobs. 


D. Solution Convergence 

Convergence of the solutions was monitored by generating convergence plots such as that seen in the Figure 5. Unix 
scripts were used to compile all of the convergence information contained in the USM3D output files into a Tecplot 
format. Each plot showed the logarithm of the residual over each iteration and the changes in the aerodynamic co- 
efficients. The criteria for a fully converged solution was for each plot to show a leveling off of the quantities under 
consideration. These plots allowed for rapid determination of whether any given solution had reached a converged 
state (Figure 5a) or if the solution had attained an unsteady solution shown by oscillatory convergence (Figure 5b). 



iteration . 

x 10 



a) Smooth convergence, a = 30deg. 


b) Oscillatory convergence, a = 48deg. 


Figure 5: Typical USM3D convergence of lift coefficient for an airfoil. NACA 4415, USM3D/SA, RE = 3 million, 

Me* = 0.2. 
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E. Post-Processing 

After solution convergence was verified, the data is processed so that total forces and moments, the spanwise distri- 
bution of forces and moments, and flow visualization may be studied. As was previously mentioned, in some cases 
an oscillatory solution develops rather than single steady-state values for the forces and moments. This has only been 
observed for 2D airfoil solutions at very high angles of attack. To handle such cases, a method had to be developed in 
order to address these oscillations. 

1. Forces and moments acting on the entire wing/configuration 

The first outputs of interest are the total force and mo- 
ment coefficients acting on a wing/configuration. A 
script was used to extract the body-axis force and mo- 
ment coefficients [ Cx CyCz Cm x C My Cm z L defined par- 
allel and perpendicular to the body coordinate system, 
and stability-axis coefficients [Cl Cd ], defined parallel 
and perpendicular to the free- stream velocity. Figure 6 
shows an example of lift coefficient results for a 2D air- 
foil and a 3D wing which used the same airfoil cross 
section. 


2. Span-wise distribution of forces on a surface 

The other output of interest is the span- wise distribution 
of force coefficients, particularly the span- wise lift coef- 
ficient. The PREDISC utility [22], was used for extract- 
ing this information. PREDISC simultaneously loads 
the grid files containing the surface grid and a converted 
TetrUSS solution file containing only surface data. Data Fi § ure 6: Comparison of lift coefficients for 2D airfoil 
extraction planes can be arbitrarily defined (see Figure (NACA 4415 > and AR = 12 rectangular wing at RE = 3 
7) and PREDISC will output surface pressure and skin million, 
friction coefficients, C p and Cf respectively, along the 

surface discretized according to a fine mesh of x/c, y/c locations. The surface pressure and skin friction coefficients 
are integrated to approximate body-axis force coefficients. An example lift coefficient distribution obtained by inte- 
grating the C p and Cf values at each extraction plane is shown in Figure 7. 





Figure 7: Screen shot of PREDISC code (left) showing defined data extraction planes and a plot of the extracted 
spanwise C/ distribution (right) on a rectangular wing AR =12. 
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IV. Results 


The results shown in this section represent the current level of development of the CFD database. The results from 
analysis of 2D airfoil CFD solutions for cambered, symmetric and thin airfoils will be presented along with 3D wing 
solutions with rectangular, tapered and swept planforms. 

A. 2D Airfoil Results 

Three different airfoils have been studied and added to the CFD aero- 
dynamic database to date. Additional airfoils will be added as needed 
for further development of the database and as required for validation 
of the low-order post-stall method discussed previously. The airfoils 
chosen exhibit different stall and post-stall behavior and offer insight 
on how certain airfoil geometries will tend to behave at high angles 
of attack. Results for a symmetric airfoil (NACA 0012), a cambered 
airfoil (NACA 4415) and a very thin airfoil (NACA 63006) are shown 
in this section. The geometries of the airfoils are shown in Figure 8. 

One interesting phenomenon that was encountered while devel- 
oping the airfoil database, was the tendency for the airfoil solutions 
to exhibit oscillatory behavior in the force and moment convergence 
histories at angles of attack of approximately 40 degrees and above. 

The cause of these oscillations was determined using flow visualiza- 
tion which revealed periodic vortex shedding from the upper surface 
of the airfoil. A process to average the oscillatory behavior and de- 
termine peak to peak amplitudes was established. A post-processing 
MATLAB script was developed to read all of the force and moment 
history files for an airfoil and detect the oscillatory behavior. For any angle of attack that displayed this behavior, the 
script identified two complete cycles at the end of the convergence history and determined the mean value along with 
the peak to peak amplitude. The plots in Figure 9 illustrate the approach used for the processing of the raw CFD data 
with this code. 



Figure 8: Comparision of airfoil geometry for 
the three cases studied. 


o 


raw cfd data 



1.5 
iteration 

2 complete cycles identified 


x 10 



Figure 9: Illustration of approach used for averaging oscillator airfoil CFD convergence history. 


8 of 18 


American Institute of Aeronautics and Astronautics 


1. Reynolds Number Effects Through Post-Stall 

The general effects of Reynolds number on the aerodynamic quantities, specifically C^ max , are known and have been 
observed with the current work as well. However, it is important to extend this to the post- stall region. Figure 10 shows 
a comparison of the lift curves for the NACA 4415 airfoil at three different Reynolds numbers (3, 6 and 10 million). 
The increased maximum lift coefficient with increased Reynolds number is expected. In the post stall region between 
angle of attack of 40 and 70 degrees, an additional effect of the Reynolds number is seen. The region falls directly 
where the oscillatory solutions develop; the values shown in Figure 10 are the averaged values from any oscillating 
data. After this recovery region, the solutions tend to follow a similar path out to 90 degrees. 

2. Sharp vs 1 . Blunt Trailing Edge Geometries 

It has been noted by Hoerner [23], that the trailing edge shape of an airfoil has a distinguishable effect on the C/ vs. 
a curve. When comparing an airfoil with a sharp trailing edge with the same airfoil but with a blunt trailing edge, it 
is found the the maximum lift coefficient is seen to be higher for the blunt trailing edge airfoil. The plot in Figure 
11 displays results for an NACA 4415 airfoil with both a sharp and blunt trailing edge at a Reynolds number of 3 
million. The blunt TE geometry was generated by removing the final one percent of the chord of the sharp trailing 
edge geometry, thus no other alterations to the geometry are present. The expected trend is seen with the blunt trailing 
edge case having a higher maximum lift coefficient. It is interesting to see that this effect seems to continue all the 
way past stall and through to approximately a = 10 degrees, after which the lift curves coincide. The flow mechanism 
that allows for this all the way through post stall is not readily apparent from the CFD at this time, but warrants further 
investigation. 




a a 


Figure 10: Effect of Reynolds number on lift coeffi- Figure 11: Effects of trailing-edge sharpness for an 
cient for an NACA 4415 airfoil. USM3D/SA, = NACA 4415 airfoil. USM3D/SA, Re = 3 million, 

0 . 2 . = 0 . 2 . 


3. Comparison of Post-Stall Characteristics of Three Airfoil Geometries 

Airfoils that exhibit different characteristics in terms of maximum lift coefficient and stall behavior were analyzed 
for addition to the post-stall database. In Figure 12, lift curves for three airfoils are shown from 0 to 90 degrees 
angle of attack. A cambered airfoil (NACA 4415), a symmetric airfoil (NACA 0012) and a very thin airfoil (NACA 
63006) are compared in the figure. The NACA 4415 and NACA 0012 both exhibit trailing edge stall behavior. That 
is, flow separation begins at the trailing edge and progresses towards the leading edge as the angle of attack increases. 
This can be seen in the lift curves as both the airfoils have a relatively “gentle” stall. It is interesting to note that 
although the NACA 4415 airfoil has a higher C/ ?majc as compared to the NACA 0012, the two airfoils have very similar 
maximum recovery lift coefficients at approximately 50 degrees angle of attack. The NACA 63006 produced a much 
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more unusual lift curve. This airfoil is categorized as having “thin airfoil” stall characteristics. Due to the presence of 
a severe adverse pressure gradient around the leading edge at even low angles of attack, the flow over the upper surface 
completely separates. This sudden onset of flow separation can been seen by the sudden leveling of the lift coefficient 
at around 10 degrees. After this initial “stall”, the airfoil recovers past the initial stall lift coefficient and reaches nearly 
the same maximum recovery C\ at a = 50 degrees as the NACA 4415 and NACA 0012. The differences seen between 
these three airfoils in the 40 to 60 degree range can mainly be attributed to thickness, camber and leading-edge radius 
effects. 

An enlarged plot of the NACA63006 in the region of incipient stall is shown in Figure 13 with comparison to 
experimental data obtained from Abbott and von Doenhoff [1] and also noted by Hoerner [23]. It should be noted 
that since the CFD solutions were fully turbulent (no transition model was used) no evidence or effects of a laminar 
separation bubble near the leading edge were modeled, a phenomenon noted by Leishman and in other studies of 
similar airfoils [24]. 




Figure 12: Comparison of C/ vs. a for three airfoil Figure 13: Comparison of C/ vs. a from CFD and 
types. USM3D/SA, Re = 3 million, M oo = 0.2 experiment for the NACA 63006 thin airfoil. 


Figures 14 and 15 show Mach number contour plots along with surface C p distributions for the NACA 4415 and 
NACA 63006 airfoils at different stages along their respective lift curves. The NACA 4415 flow visualization in Fig. 
14 at a = 17 degrees. (C^ max ) shows the presence of trailing edge stall propagating forward along the chord. This 
is also evident in the C p distribution. As the angle of attack is extended above stall and into post-stall, a large wake 
region is seen to develop accompanied by a loss of much of the suction peak. At a = 50 deg. a snapshot of the Mach 
contours and C p distribution show the unsteady behavior that exists. At a = 88 degrees the solution shows massively 
separated flow typical of bluff bodies. A small suction peak is still generated around the leading edge resulting in a 
small amount of lift generation. 


10 of 18 


American Institute of Aeronautics and Astronautics 





a) a = 17deg. 


b) a = 39 deg. 




d) a = 88 deg. 


Figure 14: Mach contours and C p distribution for a NACA 4415 airfoil. USM3D/SA, Re = 3 million, = 0.2. 


The flow visualization of the NACA 63006 CFD solutions in Fig. 15 provides an interesting view of thin-airfoil 
stall behavior. It can be seen from Figure 15 that, at a = 6 degrees, there is a significant suction peak followed by a 
steep adverse gradient. The leading edge at this angle of attack shows a small area of separated flow. Looking back 
at Figure 13, a slight change in the slope of the lift curve can be identified in this region. At 10 degrees angle of 
attack, the flow on the upper surface is completely separated; this corresponds with the maximum lift coefficient prior 
to the airfoil entering into the recovery region in post-stall. The C p distribution at 24 degrees angle of attack shows a 
larger internal area signaling an increase in the normal force generated, the majority of which still falls along the lift 
direction. The airfoil reaches another “C/, majc ” at approximately 50 degrees and also shows unsteady characteristics. 
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a) a = 6 deg. b) a = 10 deg. 



Figure 15: Mach contours and C p distribution on a NACA 63006 airfoil. USM3D/SA, Re = 3 million, M oo = 0.2. 


B. 3D Wing Results 

Results of various 3D wing geometries are presented in the following subsection. The plots in Figures 16 and 17 dis- 
play the lift curves for multiple wing planform geometries. Only a subset of the geometries have been run completely 
to 90 degrees due to the extensive run time of the solutions and the limited availability of computing nodes. The 
effects of both taper and sweep have been studied. The database will eventually include wing geometries that have a 
combination of both sweep and taper, similar to civilian transport aircraft. This data will be of interest to continuing 
efforts aimed at improving safety of transport aircraft. Detailed spanwise load distributions and surface streamline flow 
visualizations are presented in this paper for three wings of AR =12 with rectangular, tapered and swept planforms. 

1. Rectangular Wing 

The results presented are for a rectangular wing with an aspect ratio of 12 and an NACA 4415 airfoil cross section. 
The solutions for this case were produced for a Reynolds number of 3 million based on the chord and a Mach number 
of 0.2 . The plots in Figure 18 show a snapshot of the spanwise C/ distributions calculated from the time-accurate 
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Figure 16: Effect of taper on unswept wings: Cl vs. a. 
USM3D/SA, Re = 3 million, = 0.2. 


Figure 17: Effect of sweep on constant-chord wings: 
Cl vs. a. USM3D/SA, Re = 3 million, M «, = 0.2 


CFD solutions for angles of attack near stall and into post-stall. It is seen that at an angle of attack of 18 degrees, 
a sawtooth pattern in the C/ distribution is present. The extent of this sawtooth pattern seems to grow at 22 degrees 
and 28 degrees. Correlating these load distributions with flow visualization at the same angles of attack enlightened 
the reason for the sawtooth patterns. Through the use of surface streamlines it can be seen in Figure 19 that as the 
angle of attack increases, reversed flow is seen aft of the separation line and shows the presences of multiple stall cells 
forming along the semi- span of the wing. This causes certain sections along the wing to have more attached flow 
than others, generating the oscillations in the the local lift coefficients going from the root to the tip. This stall cell 
formation eventually dissipates as the flow over the upper surface becomes fully separated and the region of reversed 
flow reaches the leading edge of the wing as can be seen at 28 degrees. The streamlines in Figure 19 seem to suggest 
that even in high angle of attack situations the surface flow is still relatively in the chord- wise direction for the majority 
of the semi- span. There are some variations near the borders of stall cells and near the wing tip as would be expected 
due to the influence of the tip vortex. 




a) Local C/ distribution at various angles of attack b) Corresponding points on Cl vs. a curve 


Figure 18: Comparison of local C/ distribution at various angles of attack for a rectangular wing AR =12. USM3D/SA, 
Re = 3 million, = 0.2. 
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a) a = 14 deg. 


b) a = 18 deg. 
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c) a = 22 deg. 


d) a = 28 deg. 


Figure 19: C p contours and surface streamlines on a semi-span rectangular wing, AR =12, symmetry plane indicated 
by black line. USM3D/SA, Re = 3 million, = 0.2. 


2. Tapered Wing 

Results for a tapered wing with a taper ratio of 0.5 and aspect ratio of 12 are presented in Figures 20 and 21. This 
wing also has the same NACA 4415 airfoil cross section and the solutions were generated for a Mach number of 0.2 
and Reynolds number of 3 million based on the mean aerodynamic chord. The C/ distributions in Figure 20 show that 
the highest C/ is seen near the middle of the semi- span at 14 degrees angle of attack. 




a) Local C/ distribution at various angles of attack b) Corresponding points on Cl vs. a curve 


Figure 20: Comparison of local C/ distribution at varying angles of attack for a tapered wing, A = 0.5 and AR = 12. 
USM3D/SA, Re = 3 million, = 0.2. 

Such increases in local span C/ can be an indicator of impending stall in that region. As anticipated, a stall cell develops 
near this high C/ area of the semi-span at 18 degrees, which is denoted by the severe drop in the local lift coefficients. 
Looking at the surface streamlines in Figure 21 shows the existence of the stall cell as predicted by CFD. Similar to the 
rectangular case, these stall cells disappear once the angle of attack is high enough that the flow on the upper surface 
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is entirely separated. Chord-wise reversed flow is still seen along the majority of the semi-span at highly separated 
conditions. 



r 


4.5 -4 -3.5 -3 -2.5 


0 0.5 



a) a = 14 deg. 


b) a = 18 deg. 



CP: -4.5 -4 -3.5 -3 -2.5 -2 -1.5 -1 -0.5 0 0.5 



c) a = 30 deg. 


d) a = 42 deg. 


Figure 21: C p contours and surface streamlines on a tapered wing, A = 0.5, AR = 12, symmetry plane indicated by 
black line. USM3D/SA, Re = 3 million, = 0.2. 


3. 30 Degree Swept Wing 

Results are presented for a swept wing case with 30 degrees of sweep, an aspect ratio of 12 and an NACA 4415 
cross section parallel to the plane of symmetry. The solutions were generated at a Mach number of 0.2 and Reynolds 
number of 3 million. Span wise C/ distributions are shown in Figure 22. Stall progression from near the tip of the wing 
towards the inboard sections of the wing is evident and expected. As the angle of attack increases to 30 degrees the 
inboard local lift coefficients (calculated here based on a streamwise section, not chord- wise) are well above observed 
maximum lift coefficient seen for the NACA 4415 2D airfoil CFD solutions (approx 1.75). This phenomenon however 
has been observed in past studies of swept wings. Hunton and James [25] as well as Harper and Maki [10] note this 
same behavior, which is most pronounced at inboard sections near the root, but can occur along most of the semi- 
span. They describe this as an effect of a ’’natural boundary-layer control” that delays stall on the inboard section of 
a swept wing. It should also be noted that the percent increase in the local sectional C^ max goes up with the sweep 
angle. Evidence of this delayed stall near the root of the 30 degree swept wing can be seen in Figure 23. The inboard 
section does not show significant flow separation even at 30 degrees angle of attack, while the rest of the wing is fully 
separated. 
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a) Local C/ distribution at various angles of attack b) Corresponding points on Cl vs. a curve 

Figure 22: Comparison of local C/ distribution at various angles of attack on a swept wing, A = 30 deg., AR = 12. 
USM3D/SA, Re = 3 million, = 0.2. 



c) a = 20 deg. 


d) a = 30 deg. 


Figure 23: C p contours and surface streamlines on a swept wing, A = 30 deg., AR = 12, symmetry plane indicated by 
black line. USM3D/SA, Re = 3 million, = 0.2. 


V. Conclusions 

This paper presents on-going research related to the development of an aerodynamic database of airfoils and wings 
through post-stall angles of attack. Such data has potential for use in modeling post-stall flight dynamics of fixed- 
wing and rotary-wing aircraft, and for prediction of wind turbine performance. Except for a few sources, however, 
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there is a dearth of post-stall aerodynamic data from experimental or computational studies. One objective of the 
current work is in filling this gap in the knowledge. Further motivation behind developing this database is that it is 
expected to be of use in the validation and refinement of a low-order post-stall aerodynamics prediction method being 
developed for wings and aircraft configurations. The current approach to the development of this database is to use 
the NASA TetrUSS CFD package to analyze the geometries using Reynolds Averaged Navier Stokes equations. The 
CFD analyses enables flow solutions to be generated for various geometries (both 2D and 3D) in a more rapid fashion 
than would be possible with experimental work, besides providing detailed span wise lift coefficients, and separation 
patterns that are not easily available from experimental studies. 

In this paper, results are first presented for three airfoils that exhibit different stall and post- stall behaviors. The 
results for these and other airfoils will be used as input data for the low-order post-stall prediction method for finite 
wings. Another objective in obtaining these results is that studying the post-stall behavior of these airfoils may lead 
to the development of a rapid method of generating airfoil force and moment curves from 0 to 90 degrees angle of 
attack using a combination of results from XFOIL and other simplified models such as the flat plate theory for very 
high angles of attack. Next the paper presents results for finite wings to illustrate the effects of taper ratio and sweep 
angles on stall and post-stall behavior. The results for the rectangular and tapered unswept wings show that the flow 
along wing sections is nominally two-dimensional even at post-stall conditions. The nominally two-dimensional flow 
provides confidence in the use of sectional data in modeling post-stall aerodynamics of finite wings. However, there 
is noticeable variation in the shape of the upper-surface flow separation line, resulting in saw tooth oscillations in 
the span wise lift-coefficient distributions. It is not clear if these saw tooth oscillations have any correspondence with 
similar oscillations seen in the results predicted by the low-order post-stall method. 

In contrast to the results for unswept wings, the results for swept wings are seen to be highly three dimensional, 
as expected. At and beyond stall, there is significant spanwise flow on the upper surface. As a result of this spanwise 
flow (resulting from spanwise pressure gradients) and the higher lift coefficients on the outboard portions of the 
wing in pre- stall conditions, the outboard portions of the swept wings stall first while the root portion maintains 
lift coefficients much higher than the maximum lift coefficient in two-dimensional flow. The consequences of such 
behavior for aircraft stall characteristics, namely tip stall with associated rolling moment and pitch-up moment at stall, 
are well known. In the context of the low-order modeling, the significant three-dimensional flow on swept wings at 
and beyond stall poses serious stumbling blocks. It remains to be seen how these effects can be captured correctly 
and efficiently in a low-order aerodynamic model that can be used in real-time flight dynamics simulation. The 
benefit of such CFD studies, however, is that the detailed results do provide the very type of sought-after information 
for developing phenomenological augmentation of low-order approaches even when the flow is not entirely two- 
dimensional. In follow-on work, continued expansion of the database will progress with the addition of more airfoil 
and wing geometries. Studies of wing-tail geometries and effects of angular velocities are also planned. 
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